///@file Matrix.h
///@brief Declaration of the class representing a matrix
///@author Arnaud Duval
///@version 0.0
///@date 2010/02/25


#ifndef MATRIX_H
#define MATRIX_H

#include <iostream>
#include "Vector.h"

///@class Matrix
///@brief Class defining a Matrix
/// A matrix is an array of Vectors
class Matrix
{
	private:
		
	protected:
		///@brief Number of columns and rows
		unsigned int nbRows, nbColumns;
		///@brief Array of rows
		Vector *arows;
	public:
		///@brief Default constructor
		Matrix();
		///@brief Constructor
		Matrix(unsigned int, unsigned int, double val=0.);
		///@brief Copy constructor
		Matrix(const Matrix &);
		///@brief Copy constructor with temporary object
		Matrix(Matrix*);
		///@brief Constructor of a special square matrix
		Matrix(unsigned int, char);
		///@brief Destructor
		~Matrix();
		///@brief Overload of () operator
		double& operator()(unsigned int, unsigned int);
		///@brief Overload of [] operator : return a row
		Vector& operator[](unsigned int) const;
		///@brief Overload of = operator
		Matrix& operator = (const Matrix&);
		///@brief Overload of = operator for temporary objects
		Matrix& operator = (Matrix*);
		///@brief Copies Matrix without size check (use with caution)
		Matrix& operator &= (const Matrix&);
		///@brief Overload of += operator
		void operator += (const Matrix&);
		///@brief Overload of -= operator
		void operator -= (const Matrix&);
		///@brief Return number of rows
		unsigned int NbRows() const;
		///@brief Return number of columns
		unsigned int NbColumns() const;
		///@brief Extract a row from the matrix
		Vector ExtractRow(unsigned int);
		///@brief Extract a column from the matrix
		Vector ExtractColumn(unsigned int);
		///@brief Extract diagonal from the matrix
		Vector ExtractDiag();
		///@brief Swap two rows
		void swapRow(unsigned int, unsigned int);
		///@brief Swap two columns
		void swapColumn(unsigned int, unsigned int);
		///@brief Invert matrix
		Matrix Invert();
		///@brief Invert matrix with Gauss-Jordan method
		Matrix InvertGaussJordan();
		///@brief Invert matrix with Cholesky method
		Matrix InvertCholesky();
		///@brief Invert matrix with analytic method (limited to 2x2 and 3x3 matrices)
		Matrix InvertAnalytic();
		///@brief Transpose matrix
		Matrix Transpose();
		///@brief Compute matrix determinant
		double Det();
		///@brief Compute Frobenius norm
		double FrobeniusNorm();
		///@brief Add a submatrix into current matrix
		void AddSubMatrix(Matrix, unsigned int, unsigned int);
		///@brief Compute QR decomposition
		void QRDecomposition(Matrix &, Matrix &);
		///@brief Compute Lanczos algorithm for eigenvalues extraction
		void Lanczos(Matrix &, Matrix &);
		///@brief Check if Matrix is symmetric, assuming a tolerance
		bool IsSymmetric(double);
        ///@brief Set all components to zero
        void Nullify();

		///@brief Return true if Matrix exists
		friend bool operator ! (const Matrix&);
		///@brief Overload of + operator
		friend Matrix operator+ (const Matrix&, const Matrix&);
		///@brief Overload of - operator
		friend Matrix operator- (const Matrix&, const Matrix&);
		///@brief Overload of * operator : cross produt of two vectors A_ij = B_i * C_j
		friend Matrix operator * (const Vector&, const Vector&);
		///@brief Overload of | operator : contracted product : A_i = B_ij * C_j
		friend Vector operator| (const Matrix&, const Vector&);
		///@brief Overload of | operator : contracted product : A_i = C_j * B_ij
		friend Vector operator| (const Vector&, const Matrix&);
		///@brief Overload of || operator : twice contracted product A = B_ij * C_ij
		friend double operator|| (const Matrix&, const Matrix&);
		///@brief Overload of | operator : contracted product : A_ij = B_ik * C_kj
		friend Matrix operator | (const Matrix&, const Matrix&);
		///@brief Overload of * operator : product scalar * matrix
		friend Matrix operator * (const double&, const Matrix &);
		///@brief Overload of * operator : product matrix * scalar
		friend Matrix operator * (const Matrix&, const double &);
		///@brief Overload of / operator : division Matrix / scalar
		friend Matrix operator / (const Matrix&, const double&);
		///@brief Overload of << operator
		friend std::ostream& operator<<(std::ostream &, const Matrix&);
};

#endif

